Olfactory sampling volume for pheromone capture by wing fanning of silkworm moth: a simulation-based study

Odours used by insects for foraging and mating are carried by the air. Insects induce airflows around them by flapping their wings, and the distribution of these airflows may strongly influence odour source localisation. The flightless silkworm moth, Bombyx mori, has been a prominent insect model for olfactory research. However, although there have been numerous studies on antenna morphology and its fluid dynamics, neurophysiology, and localisation algorithms, the airflow manipulation of the B. mori by fanning has not been thoroughly investigated. In this study, we performed computational fluid dynamics (CFD) analyses of flapping B. mori to analyse this mechanism in depth. A three-dimensional simulation using reconstructed wing kinematics was used to investigate the effects of B. mori fanning on locomotion and pheromone capture. The fanning of the B. mori was found to generate an aerodynamic force on the scale of its weight through an aerodynamic mechanism similar to that of flying insects. Our simulations further indicate that the B. mori guides particles from its anterior direction within the ~ 60° horizontally by wing fanning. Hence, if it detects pheromones during fanning, the pheromone can be concluded to originate from the direction the head is pointing. The anisotropy in the sampling volume enables the B. mori to orient to the pheromone plume direction. These results provide new insights into insect behaviour and offer design guidelines for robots for odour source localisation.

Supplementary material for "Olfactory sampling volume for pheromone capture by wing fanning of silkworm moth: a simulation-based study"

S1. Kinematic reconstruction
The three-dimensional reconstruction of landmarks is described in Section 2.2.The reconstruction of wing outlines and feathering angles is described in detail in this section.In each frame, three points on the leading edge and five points on the trailing edge were manually digitised to extract the wing outline from each image.The distribution of the feathering angle along the spanwise axis of the wing was adjusted by fitting the wing root and tip of the wing morphology model to the measured threedimensional coordinates.The spanwise distribution of the feathering angle α is assumed as follows: where αb and αt are the feathering angles at the wing root and wing twist, respectively; r is the distance from the wing root; and R is the wing length.Using αb, αt, and n as variables, the distances between the model outline and the leading and trailing edges in each image were minimised using the optimisation toolbox in MATLAB (fmincon).In all frames, the average minimum distance between the model outline and the outlines in the image was approximately 3.2 pixels, indicating a sufficiently good agreement.Examples of fitting are shown in Supplementary Figure 1.Equation (S1) is suitable for the above minimisation because it enables the morphology to be expressed with fewer variables; however, n does not vary smoothly over time.Despite a good agreement being obtained near the wing tips, the relative movement of the fore and hind wings was not considered, resulting in large deviations in the outlines near the wing roots (see Supplementary Figure 1c, d).To ensure smooth variation in the parameters and enable interpolation in time, the distribution of the feathering angles was fitted again with a polynomial at each time step as follows: A feathering angle ranging from 40 % to 90 % of the wing length was used for the distribution of feathering angles defined by Equation (S1).For subsequent numerical analyses, the positional angle, elevation angle, αb, αsq, α1, α2 of the feathering angles were interpolated using a fifth-order Fourier series.The measured flapping angles and angles of attack for all individuals (k1-k3) are shown in Supplementary figure 2.

S2. Validation and verification
The simulation results were compared with the measured airflow velocities to validate the reliability of the simulations.A hot-wire anemometer (Smart CTA 7250, Kanomax Japan Inc.) was placed near the antennae of the tethered B. mori to measure the airflow induced by the flapping wings.The hot-wire anemometer was calibrated to measure air velocity (0-1 m s -1 ) with high accuracy using the voltage at the probe and measurements from a calibrated anemometer (Climomaster, Kanomax Japan Inc.) obtained in a small suction-type wind tunnel with a PC fan (Supplementary Figure 3a, b).A hot-wire anemometer and two high-speed cameras, calibrated as described in Section 2.2, were synchronised to simultaneously capture images of the flapping wings.A total of 10-57 wingbeats were recorded.The three-dimensional coordinates of the head and tip of the abdomen, wing tips, wing roots, and the probe of the hot-wire anemometer were extracted from the high-speed images to measure the position of the hot-wire anemometer relative to the body (Supplementary Figure 3c).The flow velocity in each flapping phase was calculated based on the time series of wing angles calculated from the coordinates of the wing tips.The measured speeds were filtered using a Butterworth filter with a cutoff frequency of 500 Hz.The wind speed of the fan, which was estimated using a hot-wire anemometer, was approximately 0.082 m s -1 .Subsequently, the estimated speed was subtracted from the measured speed.
Therefore, the position of the hot-wire anemometer and flow velocity were normalised by the wing length and mean wingtip velocity (2ΦRf, where Φ, R, and f are the wingbeat amplitude, wing length, and wingbeat frequency, respectively).As shown in Supplementary Figure 3d, the measured flow velocities and simulation results are in very good agreement in terms of phase and amplitude; therefore, the simulation results are reliable.
The self-consistency of the grid was tested by comparing the results from the finer mesh and time steps (finer) used in this study with those from coarser meshes and time steps (coarse) and finer meshes and time steps (finest) (Supplementary Table 2).The total aerodynamic forces and powers and the flow velocity (60 th wingbeat cycle) from each mesh are shown in Supplementary Figure 4.The flow velocity was sampled in a vertical plane 40 mm from the head of the B. mori.Despite the large differences between the coarse mesh and other two cases, the finer and finest cases showed reasonable agreement.Therefore, fine grids and timesteps were used in this study.The time-series of the flow velocity at a point 40 mm from the head of the silkworm moth is shown in Supplementary Figure 4c as an example of the time history of the far-field flow.As the flow velocity is well converged after approximately 60 cycles, the repetition of the flow field at the 60th cycle was used to calculate the particle trajectories.

Supplementary Figure 4
Verification of the computational grid and flow field convergence.(a) Time series of the resultant aerodynamic force and power; and (b) flow velocity simulated by coarse, finer, and finest meshes and time steps.The flow velocity in (b) was sampled at the vertical plane placed 40 mm from the head of the silkworm moth model.In the central diagram in (b), the a silkworm moth model and the position of the flow velocity are indicated by red dots, and the corresponding figure of time series of the flow velocity is connected to the red dots by black lines.(c) Time-series of flow velocity at the furthest point from the model (the blue circle in (b)) among the sampling points in (b).Supplementary Figure 5 Time series of (a) aerodynamic force and (b) power for each individual (k1-k3).The weights of each individual are shown by horizontal lines for comparison.Horizontal, side-slip, and vertical forces represent the aerodynamic force along the x, y and z axes in Figure 1, respectively.The side-slip force is cancelled when the left and right wings are combined; therefore, the side-slip force for one wing is shown in the figure.Supplementary Figure 7 Time series of the violin plot of the distance (top), horizontal (middle), and vertical angles (bottom) of the particles with respect to the middle of the antenna.The arrival time of the particles at the antenna is set to zero.